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ABSTRACT 

We analyse some properties of circumplanetary discs. Flow through such discs may 
provide most of the mass to gas giant planets, and such discs are likely sites for the 
formation of regular satellites. We model these discs as accretion discs subject to the 
tidal forces of the central star. The tidal torques from the star remove the disc angular 
momentum near the disc outer edge and permit the accreting disc gas to lose angular 
momentum at the rate appropriate for steady accretion. Circumplanetary discs are 
truncated near the radius where periodic ballistic orbits cross, where tidal forces on 
the disc are strong. This radius occurs at approximately 0.4 rn for the planet Hill radius 
7'H. During the T Tauri stage of disc accretion, the disc is fairly thick with aspect ratio 
H/r > 0.2 and the disc edge tapering occurs over a radial scale H ^ O.lrn. The 
disc fluid equations can be rescaled in the Hill approximation to a form similar to 
the flow equations for a disc in a binary star system with a mass ratio of unity. For 
a circular or slightly eccentric orbit planet, no significant resonances lie within the 
main body of the disc. Tidally driven waves involving resonances nonetheless play an 
important role in truncating the disc, especially when it is fairly thick. We model the 
disc structure using one dimensional time-dependent and steady-state models and also 
two dimensional SPH simulations. The circumplanetary disc structure depends on the 
variation of the disc turbulent viscosity with radius and is insensitive to the angular 
distribution of the accreting gas. Dead zones may occur within the circumplanetary 
disc and result in density structures. If the disc is turbulent throughout, the predicted 
disc structure near the location of the regular Jovian and Saturnian satellites is smooth 
with no obvious feature that would favor formation at their current locations. It may 
be possible that substructure, such as due to variations in the disc turbulence, could 
lead to the trapping of migrating satellites. 

Key words: accretion, accretion discs - planets and satellites: formation - planetary 
systems - planet-disk interactions 



1 INTRODUCTION 



In the core accretion model of planet formation, at t he earliest stages of giant planet formation, the planet is fully embed- 
ded in the gaseous disc that orb its the central star ( Mizunol 1980l : Pollack et al. 19961 : Hubickvi. Bodenheimer fc Lissauer 



20051 : iPapaloizou fc Nelsonlboosi ). At later stages, during run-away-gas accretion, tidal forces due to the planet open a gap 
in the disc, typically when the planet reaches a mass of order Neptune's mass (jPin fc Papaloizoul 1 19861 : iBate et al.l bood : 
P'Angelo. Henning fc Klevlbood') . The gap opening does not necessarily imply that the gas flow on to the planet has ceased 
jArtymowicz fc Lubowlll996l 'l. The gap structure is determined by the properties of dynamical flow within it. The relatively 
low density of gas in the gap results from the relatively high flow radial velocities towards the star as a consequence of the 
disturbance caused by the planet. Outside the gap region, the gas radial drift occurs on the much slower viscous timescale, 
resulting in much higher densities than occur in the gap. Some studies suggest that most of the inflowing gas just outside 
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the Rap is accreted by a planet whose mass is comparable to or less than Jupite r's for typical disc parameters (jBrvden et al 



the gap IS accreted by a planet whose mass is comparable to or less than Jupite : 
I999I : lKlevlll999l : iLubow. Seibert fc Artvmowiczlll999l : iLubow fc D'AngelollioO^ ) 



During this gap phase, the planetary radius is typically much smaller than its Hill (tidal) radius. Cir cumplanetary discs 
may form as gas flowing through the g ap on to a planet carries some angular momentum about the planet ( Lubow et al. 19991 : 



Bate et al. 20031 : D'Angelo et al. 2002). As a result, the inflowing gas has too much angular momentum to directly strike the 



planet and instead forms a disc about it. Most of the mass of Jupiter may have been acquired in the gap stage through gas 
flow that involves a circumplanetary disc. 

There are at least two motivations for studying circumplanetary discs. The observational detection of a disc orbiting 
an extra-solar gas giant planet would provide important evidence about the process of planet formation. The determination 
of their expected properties provides constraints for their detection. A second motivation is to better understand satellite 
formation, since circumplanetary discs are likely sites of satellite formation. Satellites that orbit Jupiter and Saturn are 
classified into two groups. The regular satellites have low eccentricity and low orbital inclination to the equatorial plane of 
their planet. They are thou ght to have formed in the circumplanetary disc ( Lunine fc Stevenson 19821 : Canup fc Ward 2002 : 



Mosqueira fc Estrada|[20o3 ) . The second type, the irregular satellites, can have a high eccentricity a nd high orbital inclinations. 



They can orbit either progradely or retrograd ely with respect to the planet sp in dCrav et al.l|2003! ). The circumplanetary disc 
may capture these satellites though gas drag JPoUack Burns fc Taubeij[l97alCuk fc Burnsl l l2004f ) , although the disc does not 
play such a role in the Nice model (e.g. iGomes et al.ll2005l ). The circumplanetary accretion disc plays a vital role in satellite 
formation. The regular satellite systems of Jupiter and Saturn extend over a small fraction of their respecti ve Hill radii, to 



less than 0.06rH. This fact has mo tivated a dynamical explanation in terms of circumplanetary disc sizes (e.g. lCanup fc Ward 
2OO2I: iMosaueira fc Estradalliooil ) . We discuss this in more detail in Section 8. 



Previous studies have considered the role of the angular momentum (about the planet) of the inflowing gas in determining 
the disc structure, in particular the disc outer radius. This picture has some intuit ive appeal, since n o disc would form, i.e., 
the disc radius would be zero, if the inflowing gas has zero angular momentum. Quillen fc Trilling ( 19981 ) suggested that 
after the planet has opened a gap, the circumplanetary disc extends to r « rn/S . This result is based on the initial angular 
momentum of gas whose flow relative to the planet is slow near the planet's Hil l sphere and sub s equen tly a ccelerates inward 
towards the planet while conserving its angular momentum about the planet. ICanup fc WardI (I2OO2I ) and IWard fc Canup 
I 2010|)_ considered ca ses where the angular momen tum of the inflowing gas is considerably smaller. 



Lissauer (1995) (see also Lissauer et al. 20091 ) suggested that the circumplanetary disc structure prior to or during gap 



opening could account for the location of regular satellites. This model was also based on the concept that the disc structure 
is determined by the angular momentum of infalling gas about the planet, with an implied smaller value of disc radius of 
r = rH/48. This paper concentrates on the former case, the case of a circumplanetary disc in the presence of gap opening. 
The results of this paper may also apply to this case. We briefly discuss this point again below equation Q and in Section [S] 

This radius characteristic of the angular momentum of the accreting gas is useful for providing an estimate of the 
conditions required for circumplanetary disc formation. For a disc to form about a planet with a gap, the planetary radius 
must be smaller than this value. Violating this constraint for disc formation would require very short orbital periods, of order 
1.5 d, for a planet of Jupiter's radius. 

However, it is not clear that this radius characteristic of the angular momentum of the accreting gas is important for 
determining the disc radius or even has much influence of the disc structure. The reason has to do with the requirements of 
angular momentum loss from a steady-state di sc. The process o f angu lar momentum transport in circumplanetary discs has 
up to now been largely ig nored. The work by ICanup fc WardI tooi) and IWard fc Canup I |2O10h did consider the angular 
momentum issue, as will be discussed further in Section 6.5. IMosaueira fc Estradal (|2003l ) briefly discuss the possible role 
of resonances in truncating the disc. For the disc to accrete in a steady-state manner, gas within it must continuously lose 
angular momentum. These studies recognized that the disc would extend beyond the radius in which mass is injected in 
order to account for the requirements of angular momentum conservation. T he angular momentum issue has been previously 
encountered in the case of discs in mass exchange binary star systems (e.g. iPapaloizou fc Pringlelll977l ). The inflowing gas 
from the mass losing star is captured within the Roche lobe of the companion star. This gas has some angular momentum 
that often results in the formation of disc about the companion. The picture that has resulted is that the disc spreads outward 
by viscous torques until it reaches a radius at which tidal torques from the companion are sufficient to remove the angular 
momentum at the rate required for steady state accretion. The disc structure is insensitive to the angular momentum of the 
accreting gas. The case of a planet accreting gas can be considered to be a binary of extreme mass ratio. The accretion disc 
model in the binary star case suggests a somewhat different picture of how the accretion process operates than has been 
previously considered for circumplanetary discs. We consider the accretion disc model of circumplanetary disc s in this paper. 

S ome recent three-dimensional simulations have analysed the flow in circumplanetary discs. For example, lAvliffe fc Bate 
( 20091 ) simulated gas accretion by protoplanets in three-dimensions with gas self-gravity and radiation transfer to investigate 
the properties of circumplanetary discs. They determined the disc scale height to be H/r > 0.2 and found the dis c size t o be 
about 0.35r-H. They attributed this value of the disc radius as following from the estimate of Quillen fc Trilling 1 19981 ), for 
which the radius is determined by the angular momentum of the accreting gas. 
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Machida l|2009l ) also performed three-dimensional calculations of circumplanetary discs. For a Jupiter mass planet, he 



found a peak in the surface density at about 21 protoplanetary radii (or 0.028 re) from the protoplanet centre. According 
to these results, the peak is due to the balance between the centrifugal force and the gravity of the protoplanet. Such small 
scale features might then correspond to the formation sites of the regular s atellites around Jupiter and Saturn. However, 
this peak in the surface density was not reproduced by Avliffe fc Batel ( 20091 ) who found that the surface density decreases 
monotonically with radius from the protoplanet. 

We analyse the general properties of circumplanetary discs by the use of some simplified models, in order to understand 
the physical processes that affect disc structure. In Section [2] we estimate the scale height and temperature of the disc. In 
Section [3] we consider a ballistic disc model consisting of periodic particle orbits around a planet. In Section 0] we examine 
the circumplanetary fluid disc evolution equations in the Hill approximation and discuss their relationship to the binary star 
case. In Section[S]we consider whether resonances could lie within a circumplanetary disc. In Section[Bl we analyse the effects 
of a tidal torque on a purely viscous disc. In Section [7| we describe results of a model that includes the effects of gas pressure 
by means of some SPH simulations. Section [8] contains a discussion and the conclusions. 



2 CIRCUMPLANETARY DISC PROPERTIES 



We consider a system with a planet of mass A/p in a circular orbit around a star of mass AIs at a separation a. The Hill sphere 
is the approximate region where the planet's gravity dominates that of the star and its radius is given by 

where = Mp/{Mp + Ms) ^ 1. The sound speed in an ideal gas is 

(2) 

fimmji V iU" 

where k is the Boltzmann constant, T is the temperature, fi^i is the mean molecular weight and mu is the mass of a hydrogen 
atom. The surface temperature of a steady state accretion disc around an object of mass M at a radius r is given by 

3 GMM 



Stt 



(3) 



( PringldllQSlI ) where a is the Stefan-Boltzmann constant. 

We estimate the temperature of the accreting gas with the disc. Similar approaches were taken bv lCanup fc WardI (|2002l ). 
Consider a disc orbiting a solar mass star that is accreting gas at a rate of lO~*M0yr~^, a typical value inferred from 
observations of T Tauri stars. Equation [3] implies that at a distance from the central star of 5 AU, the disc has a temperature 
of about 26 K, ignoring effects of stellar heating. This temperature corresponds to a disc aspect ratio H /r ~ 0.04. 
The aspect ratio of the circumplanetary disc taken with respect to the planet is 



H 



= 0.3 



rH/3 



Mr., j I Mj j I 



a \ 8 
.5AU/ 



10-« Mq yr- 



(4) 



where Mj is the mass of Jupiter and r is the distance from the planet. (There is a small correction due to the vertical 
gravity of the star that we ignore.) This aspect ratio is much larger than the value for the circumstellar disc, {H/r)s — 0.04. 
We have taken the disc sound speed Cs to be given by its value at the disc surface and have ignored any increase at the 
midplane with optical depth. The relatively high value of {H/r)p is due to the higher gas temperatures and the typically 
smaller fir values in the circum planetary disc case. T he value of the disc aspect ratio from equation ^ is similar to the values 
obtained in the simulations by Avliffe fc Batel ( 20091 ). They found no evidence for circumplanetary discs in their simulations 
of low mass planets (less than 100 Earth masses). They attributed this result to the weakness of the planet's vertical gravity 
compared to pressure, as can be seen in equation (|4]). For a sufficiently low mass planet, {H/r)p can be large enough that the 
concept of a disc breaks down and the gas around the planet is essentially spherical. As is well known, the circumplanetary 
disc temperatures at these accretion rates are too high to explain the exist ence of the icy s atellit es of Jupiter and Saturn 
I Canup fc W ard 2002). However, in a slower accretion phase as proposed bv^ Canup fc WardI (|2002l ) and further investigated 
bv lBarr fc Caimp. (2008 1. the temperature would be lower and H/r would be smaller, H/r ~ 0.1. We discuss this further in 
Section m 

We compare the surface density in the circumplanetary disc, Ep, to the surface density in the circumstellar disc, Es. We 
assume that most of the mass being accreted through the circumstellar disc outside the orbit of the planet is accreted by it. 
As discussed in the Introduction, this efficient gas capture has been found in simulations of Jupiter ma ss planets that orbit 
solar mass stars. The mass transfer rate through a steady state disc is given by M oc i/E (|Pringlelll98"ll ). where the viscosity 
is parametrised with the a prescription so that 



u — aCsH — a \ — 1 



r ii 



(5) 
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I Shakura fc Sunvaev|[l97i ) where Q — ^ GM^/r. It then follows that 

Ss V-p ,2/3 Q^p \ r ) p 



where and i^p are the viscosities in the circumstellar and circumplanetary discs respectively, and Qs and Op similarly refer 
the a values in the discs. 

Consider a solar mass star and a Jupiter mass planet, = 10~^, and circumstellar and circumplanetary disc aspect ratios 
of 0.04 and 0.3 respectively, as discussed above. For r = 0.2rH and equal a values (as = Qp) in the above equation, the 
ratio of the surface densities is Es/Sp ~ 0.25. So the circumplanetary disc surface density is somewhat higher that the local 
circumstellar disc density. Circumstellar discs could under some circumstances contain dead zones, regions where the disc 
ionization is too low for the magneto-rotational instability to operate as a source of turbulence (Gammie 1996). Dead zones 
operate where the temperature is sufficiently low that thermal ionization is weak (less than about 10^ degrees) and where the 
surface densities are high enough that external sources of ionization such as cosmic rays do not penetrate far below the disc 
surface. For expected conditions in a circumstell ar disc, dead zo nes may extend from a few tenths of an AU to several AU, in 
some cases beyond Jupiter's orbital radius (e.g. Terauem 20081 ) . The model here suggests that the conditions for dead zone 
formation within a circumplanetary disc are perhaps more favorable than in the nearby circumstellar disc gas. The reason is 
that the circumstellar disc surface densities are higher than in the nearby circumstellar disc gas, while the temperatures are 
low enough to avoid sufficient ionisation in much of the circumplanetary disc. A difference in the properties of turbulence in 
the two discs violates our assumption of Qs = Qp. But, it does so in a way that may further enforce the importance of dead 
zones in the circumplanetary disc. That is, as/op would be expected to increase, leading to a decrease in Es/Sp. We do not 
pursue the possibility of dead zones further in our disc models described later and assume a simple viscous disc. 

The viscous timescale in the disc is 

. r! . . ,„.p (iS:!) Vi) "\ 

Vp a(H/r)lQ, V " / \H/rJ \rH J 

where P is the orbital period of the planet. So we expect that such discs should be viscously relaxed for Jupiter and Saturn. 
This result has the implication that the viscous disc flow can be regarded to be in a steady state. 

A crude estimate of the steady state luminosity ratio of the circumplanetary to circumstellar discs in a system with a 
Jupiter mass planet is 

Lp ^ MpMpR, 
is ~ MsMsi?p ' 

where Rp and _Rs are the planet and star radii, respectively. For Mp ~ 1 — lOMs, Mp ~ O.OOlAfs, and Rp ~ 0.1-Rs, we estimate 
Lp/ Lb to be of order a few percent. We see that circumplanetary discs are not very bright. We discuss this further in Section[8l 



3 BALLISTIC PARTICLE PERIODIC ORBITS AROUND THE PLANET 

In this section we consider ballistic particle orbits in a planet-star system in order to find their nonKeplerian angular velocity 
and where the orbits begin to cross. 



3.1 Ballistic Equations 

As a simple two-dimensional disc model, we consider ballistic particles orbiting a planet of mass Mp in the corotating frame of 
the star-planet system with a star of mass Ms. The particles lie in the star-planet orbit plane. For such a model to represent a 
low pressure (cold) steady state fluid disc, each orbit must be periodic in the corotating frame and be nonintersecting, either 
with itself or neighboring orbits, in order that the velocity be single-valued in space. In addition, the orbit must be stable. 
The equation of motion of a ballistic particle at position r for potential (j) in the corotating frame is 

r + 2np X r = -V</>, (9) 

where fJp is the angular velocity of the planet. The potential due to the rotation of the frame and the gravity of the planet 
and the star is given by 

^_ |r-rcm|^»g GMp GM, 



^ I ' 'Pi I ' ' s 

where rem, rp and rs are the position vectors of center of mass, the planet, and the star respectively. The angular velocity of 
the star-planet binary system is 

fG{Ms + Mp, 



. ■ (11) 
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We adopt units where a — Ms + Mp — Qp — 1. The equations of motion for the particle with Cartesian coordinates {x, y), 
whose origin is at the center of mass, are then 

x^2y-x+ ^["-(1-^)] 3 + + =0 (12) 

{[x-{l-fi)f + y'^)' [{x + + y^]2 

and 

y + 2x-y+ r+ 3 =0. (13) 

Following the standard procedure to obtain the Hill equations, we transform these equations to a coordinate system 
centered on the planet and rescale the coordinates. The scaling is chosen so that the new coordinates scale with the Hill 
radius, see equation We make the change of variables to rescaled radius R defined by 

fl= (r-rp_„,)/(aM*), (14) 
where rp„cm is the displacement of the planet from the center of mass. The equation of motion @ is then 

i?+2e, X -V$, (15) 

where <I> is the potential in this new frame that we determine in Section [3.21 We let X = 1 — fi + fi^x and Y = fi'Sy. Since fi 
is small, we consider here only the terms to lowest order in ^. The equations to order /iS become 

-3X+ =--2y + X = (16) 

(X2+F2)l 

and 

3- + 2^ + ^ = 0. (17) 

(X2 + Y-2)i 

In Section [3.41 we consider the higher order terms in /i that have been neglected here. 
3.2 The Potential in the Hill Approximation 

The terms in equations p6[l and p7p that are functions of X and Y are expressed as the potential gradients 

|i = _3X+ ^ (18) 

dX (X2 + F2)f 



and 

9$ Y 



dY (x2 + y2)i 

Integrating these, we find the potential in the Hill approximation to be 



(19) 



$ = r - (20) 

(X2 + F2)i 2 

In polar coordinates centered on the planet, so that X = RcosO and Y = i?sin^, the potential is given by 

$ = - -i^^cos^e* 
R 2 

= - -i?' - -i?^cos26l. (21) 

The first term is a point mass potential of the planet. The other two terms are due to the rotation of the frame and the 
gravitational effects of the star. 

In dimensional form, the potential is given by 

^2 2 fa . 3 /r\2 
<?)H = — i " " 

which is valid for r ~ 0(/i^''^a) 



= -J^pfl M ( - + ^ (-j (l + cos2e)), (22) 



3.3 Angular Velocity in the Hill Approximation 

From the r-component of the equation of motion ((9)l with = (/)h, we have that 



^ -2f7pr(f7-fip) = (23) 
r or 



We substitute the axisymmetric terms in the potential of equation (|22|) and the velocity in the 6'-direction, ue — r(f2 — f2p) 
to find 
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Figure 1. Left: Periodic orbits about a planet in the Hill approximation. The x and y coordinates are in units of the Hill radius rjj- At 
a distance of about 0.4rH from the planet, the orbits begin to cross. Right: Stability parameter (Henon 1965) a — 1 plotted against xq, 
the X coordinate in units of rjj where the orbit crosses the positive x-axis. For a — 1 negative (positive) the orbits are dynamically stable 
(unstable). The dashed line is for periodic orbits that eventually intersect each other at larger xq, as occurs in the left panel. The solid 
line tracks the stability of periodic orbits that remain symmetric about the y-axis for larger xq. Notice that for xq ~ 0.4rH, vfhere orbits 
begin to cross in the left panel, there is a bifurcation in stability properties. At larger radii, y-symmetric orbits remain nonintersecting, 
but are unstable, while the intersecting orbits are marginally stable. 



where r is the (unsealed) distance from the planet. There are corrections to this result due to nonlinear effects of nonaxisym- 
metric terms. The angular velocity about the planet, 57, is less than Keplerian because of the effects of the star. 



3.4 Crossing and Unstable Orbits 

We numerically solved equations p6[) and ((T7|) to find the tr ajectories of part icles in this potential. We found closed periodic 



particle orbits in the frame of the planet, along the lines of iPaczvnskil (|l977l ). but in the Hill approximation. We chose the 



initial position for the particle to be on the a;-axis. The initial y- velocity is determined so that the particle has zero velocity in 
the X direction when it again crosses the x-axis on the other side of the central planet. We find that the neighboring particle 
orbits first intersect for the initial position on the a;-axis of a; = 0.410rH (or X — 0.294), as shown in the left panel of Fig. [1] 
This maximum radius of nonint ersect ing orbits, rmax — 0.410rH, is considered to represent a maximum disc radius for a cold 



disc, as discussed by Paczynski l|l97it ) 



The orbits in in the left panel Fig. [T] are asymmetric (lopsided) with respect to the y-axis in that the intersections occur 
for a; < 0. On the other hand, the m — 2 Hill potential is symmetric about the y-axis. A more detailed analysis shows that 
there is another set of orbits for which the intersections occur for x > 0. In addition, there is a third set of orbits that remain 
symmetric about the y-axis and nonintersecting for r > rmax. However, this set of orbits is dynamically unstable for r > rmax. 
That is, there is a bifurcation that occurs at r = rmax, as seen in the right panel of Fig. [T] Beyond this radius, there continue to 
exist y-symmetric nonintersecting orbits, but they are unstable. There are stable orbits beyond this radius, but they intersect 
with neighboring orbits. In any case, this critical radius of rmax ~ 0.410 ru is the limiting radius for orbits that represent a 
cold, steady-state disc. 

To test the validity the Hill approximation, we also solved the equations of motion for the full gravitational potential in 
the frame of the planet. We determined the particle orbits in the same way as described above. We found that the orbits first 
intersect at a radius of 0.412 rn with a mass ratio fi = 0.01. The Hill approximation is then reasonably accurate for these 
purposes, even for relatively large mass planets. 

3.5 Higher Order terms 

By including terms of 0(^3) and 0(/i) in equations (|16p and (|17[). we obtain the dimensionless potential 

$ = - - - -i?^ - -i?^ COS 261 + —7?^ (3 cos 6* - 5 cos 361) - ^i?^ (9 + 20 cos 29 + 35 cos 46*) . (25) 
i? 4 4 8 64 

These higher order terms in fi^^'^ introduce additional angular dependencies. 
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4 DISC FLUID EQUATIONS 

The steady-state two-dimensional equations of motion for a fluid disc in the corotating frame are 

V • (Eu) = (26) 

and 

(u-V)u + 2r2p xu = -V</)-ivp + /v, (27) 

where u is the fluid velocity and /v is the force per unit disc mass that represents the effects of turbulent viscosity. Quantity 
p is the two-dimensional (vertically integrated) pressure. We consider a cylindrical coordinate system (r, 6) centered on the 
planet. The azimuthal viscous force due to shear is given by 

1 d f sdQ\ 



where S is the surface density and v is the kinematic turbulent viscosity. The potential (f> is given by equation (|10|) . 

We rescale the variables by /i^^^ as we did for the ballistic case in equation ([T4|. We rescale the velocity and radius as 
U = u / {jj,^ aQp) and R = r/{fi^^^a), respectively, and obtain 

V ■ (EU) = (29) 

and 

(U ■V)U + 2npXU = - ^VP + F^, (30) 

where 

F. = /v/(M'^'a2r!p), (31) 

V = (M'/'a)V, (32) 

P^p/{^^/'a'nl), (33) 

n^n/Qp, (34) 

and 

$ = <^/(^i^/^a^n^). (35) 

The potential, <!?, is given by equation (|25[) and is independent of /i to lowest order, the Hill approximation. The rescaled 
viscous shear force is 

d ( f hV -dO. 



where we used the a-prescription for the viscosity so that 



H\ _ 2/3 / H\ 



—J rn^^l'■'•'a\^—j R-'VtaO.p (37) 
( Shakura fc Sunvaevlll973l ). We have left the disc aspect ratio about the planet in the unsealed form H /r because it is invariant 



2 



1 R^n^, (38) 



to scaling. We have not rescaled the surface density because doing so would have no effect on the velocity field, since we neglect 
disc self-gravity. Notice that Fr.e is independent of the mass ratio /i. This is also true of all other components of the viscous 
force. The rescaled pressure is given by 

r 

which is independent of jj.. Thus, the transformed equations ([29} and (|30p are of the same form as the original equations 
(|26|l and p7p and are independent of /i in the Hill approximation. This result also holds in an obvious extension to three 
dimensions. 

The rescaled circumplanetary disc equations are then very similar to the circumstellar disc equations previously analysed 
for close binary stars with order unity mass ratio. One difference is the form of the potential, which has a single nonaxisym- 
metric azimuthal number m = 2 in the Hill approximation. However, this is the dominant tidal term in close binaries. Consider 
a disc in binary that orbits star 1 and is tidally perturbed by star 2. At a small distance r <Si a from star 1, the potential 
terms due to the star 1 and the m — 2 tidal component due to star 2 is given hry (j) — —GM\/r — SGAh cos{29)r^ /{ia"^). 
This potential is similar (within a factor of 2) to the rescaled potential of equation (|2ip with G = l,a = 1, and masses 
All ~ M2 = 1/2. In particular, the ratio of the m = 2 tidal to central potentials is the same for both cases. Therefore, the 
disc flow equations for a planet in the Hill approximation are similar to those for a binary star system with unity mass ratio, 
for a given value of H/r and a. 
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5 RESONANCES 



Resonances occur in the disc where the forcing frequency matches a natural frequency in the disc. Angular momentum can 
be transferred to the star fr om the disc by the tidal torques that are exerted at the Lindblad resonances within the disc 
i Goldreich fc Tremainelll979l '). At such resonances, rotationally modified pressure waves are launched. Torques are exerted on 



the disc at radii where the waves damp. Resonance torques could play a role in truncating the disc. In this section we consider 
whether circular, eccentric, or vertical resonances could lie within the disc. We only consider nearly circular orbits in this 
analysis. We find this approximation holds well out to the radius where orbits begin to cross fSection l3.4|l . Since we expect 
the disc to be truncated inside or at this radius, we should find all of the resonances present in the disc in this approximation. 

Ifowever, for a disc as warm as suggested during the T Tauri accretion phase, H/r ~ 0.3 (equation Q), some gas might 
extend beyond the orbit crossing radius and off-resonant forcing could play a role. Off-resonant forcing is possible because for 
m of order unity, the resonance width scales as {H /r)^^^ m. With such a large width, the resonance could overlap with the 
disc, even though t he exact resonance location does n ot lie within the main body of the disc, as has been investigated in the 



binary star case bv lSavoniie. Papaloizou fc LinI (|l994l ) 



5.1 Circular Resonances 



We apply the angular velocity of the disc in the Hill approximation given in equation (|24p to the Lindblad resonance condition 
Goldreich fc Tremaind (|l979h . The circular Lindblad resonances then occur where 



(39) 



In the Hill approximation we have only the m = 2 term of the potential. We find that the only positions inside the orbit 
crossing radius is for the m = 1 term which is at i? = 0. For a cold disc, the resonance width is very small and the tidal forcing 
is very weak near the disc center. For a warm disc, stronger resonant excitation is possible. At higher order in ia}'^^ , other 
m- values are present (equation 125)1. However, these higher order resonances also fail to lie inside the orbit crossing radius. 



5.2 Eccentric Resonances 



resonances 



If the orbit of the planet is slightly eccentric with eccentric ity e, the tidal forcing can b e decomposed into a series or rigidly 
rotating potentials at various frequencies tQ.-p for integer t ( Goldreich fc Tremainel 19791 ). The eccentric Lindblad 
then occur where 

r dn' 



2 dr 



= 0. 



(40) 



For m — 2, the lowest order resonance that lies within the orbit crossing radius has 1 = 7. The resonant torque scales as 
g2|i-m|^ where m is the azimuthal wavenumber of the tidal potential. Consequently, the torque would scale as e^°, which 
suggests that the torque would be quite weak for modest eccentricities. For such a resonanc e to be able to overcome th e 
effects of disc turbulent viscosity and truncate the disc, we roughly require a < (r/H)'^ e^" ( 
Consequently, the disc viscosity would need to be very small. 



Artvmowicz fc Lubowlll994l ) 



5.3 Vertical Resonances 

There is also a set of resonances associated with the vertical disc motions. For each m, the vertical disc resonance lies closer 
to the planet than the corresponding horizontal (coplanar) Lindblad resonance, although the torque it produces is a factor of 
{H/r)^ smaller than the corresponding horizontal resonance. The vertical resonance generates horizontally propagating waves 
which travel towards the planet at the centre of the disc. The wave generation can transport angular momentum to the orbital 
motion of the planet-star system. The vertical resonances occur where 

m\Q-n^f = {i + r)(nl + ^], (4i) 



where F is the effective adiabatic index ( Lubow 198ll 'l. We solve this equation with F = 1.4 and find the m — 2 resonance 



occurs at r/rn = 0.49, which is somewhat outside the orbit crossing radius. Resonances at higher m-values occur further away 
from the planet. The analysis in subsequent sections of this paper considers only one- and two-dimensional models that are 
not capable of finding vertical resonances. 
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6 VISCOUS DISC MODEL 



We consider the effect of the tidal perturb ation on a circumplanetary accretion disc due to the presence of the star by following 
the work of IPapaloizou fc Pringld ( 19771 ). The tidal t orque transfers angular momentum fro m the disc to the binary orbit, 
allowing material to accrete on to the central object ( Borner et al. 19731 : Lin fc Pringle 19761 ). The disc attempts to expand 
by viscous forces, but tidal torques dominate the transport process only in the outer parts of the disc. In this approach, the 
disc pressure forces are ignored. Consequently, this approach does not include the effects of off-resonant waves that could 
be present. It is more accurate for cool discs, H/r < 0.1, where the effects of waves are likely less important. Such cool 
circumplanetary discs could arise at late stages of disc evolution where the accretion rates are lower as the disc disperses. 



6.1 Linearised Velocity 

We consider linearised equations for the tidal disturbances about the circular motions around the planet and initially include 
only gravitational and centrifugal forces. The unperturbed gas orbits with velocity uo — (0, r(f7 — ilp)) in the corotating frame 
centered on the planet, where Q is given by equation (|24p . We consider perturbations to this flow by the presence of the star. 



We de note the perturbed velocity as ixi, where ui ^ uo and ui = {ur,ue). We follow the approach of Papaloi zou fc Pringle 



(19771), but adopt the corotating frame with the Hill approximation. The linearised equations of motion based on equation (|27|l 



{n ^ Q^)^ - 2nug = fr (42) 



(n-np)^+2B«, = (43) 



and 



where 



The force components are obtained from the negative gradient of potential in equation 
Since only m = 2 disturbances are involved in the Hill approximation, we represent 

Ur,ug,fR,fe oc exp(2j6l) (45) 

and obtain 



and 



(46) 



_ -Bfr + i{n - np)/9 

* 2[-(n ~ Qp)^ + HQ]' ^ ' 



6.2 Dissipation 

We now consider the additional dissipation in the disc due to tidal perturbations when we introduce a small amount of 
dissipation into the flow. The two modes of dissipation in a Newtonian fluid are the shear viscosity, v, and the bulk viscosity, 
The dissipation due to the bulk viscosity per unit radius is 

\2 



Di = EC / (V ■ uY de (48) 
Jo 

which we can evaluate using equations (|46|) and (|47() . The dissipation due to the shear viscosity is 

r2TV 

Da = 2Si. / E,jE,j dO, (49) 







Cij — l/35ij'V .u and eij is the rate of strain tensor which is given by 
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Figure 2. Left: The dissipation in the disc. The soUd hne is Do/{uT,), the dashed line is Z)i/(fE) and the dotted line is D^Ky'i^)- Right: 
The torque acting on the disc per unit surface density, per unit viscosity. The solid line is for the viscous torque in the disc, To/{!^S), 
the dotted line is the torque Ti/(fS) and the dashed line is T2/(i/S). For both plots the vertical dotted lines show where the particle 
orbits cross. 



6.3 Torque on the disc 

From the dissipation we can find the torque on the disc due to the tidal perturbation. The rate of working per unit surface 
area of the disc by a torque, T, is 



dT 

D = n^dr 



n dT 



dr ATvrdr 47rr dr 



1 



djTQ,) 
dr 



(51) 



where Q' — d^l/dr. The first term is the rate of convection of energy over the whole disc. Its value depends only on the 
boundary conditions. We find the torque on the disc to be 



(52) 



where i = 0, 1, 2 dPrank. King fc Rainel2002h . The torque To is defined as the torque on an accretion disc without a companion. 
This has a viscous torque of 



To = 2nriyT.r^n' 



(53) 



and so the dissipation is 



Do = -r'{n'[ 



(54) 



In the left panel of Fig.[2]we plot th e three scaled dissipations, Do /(i^S), Di/(("E) and D2/(i^E) as functions of the radius 
in the disc. This is similar to Fig. 1 in iPapaloizou fc Pringld (|1977l ). but we use the Hill approximation. We see that the 
magnitude of the dissipation from the internal viscous torques, Do/{uYl), and the tidal dissipation, D2/(i^E), are equal very 
close to, but just inside, the radius in the disc where the particle orbits begin to cross. In the right panel of Fig. [2] we plot the 
scaled torques on the disc. The tidal torque in the disc starts to dominate the viscous torque just inside of the radius where 
the particle orbits cross. 

However, outside of the radius where the particle orbits cross, pressure and nonlinear effects cannot be neglected. There- 
fore, this linear solution is not valid beyond the orbit crossing radius. 
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6.4 One-Dimensional Simulations 

We use a one-dimensional (radius only) model of an accretion disc subject to a tidal torque in order to determine the surface 
density evolution. The governing equation for the surface density of a flat accretion disc centered on the planet is 



9E _ ld_ 

dt r dr 



+ S{r) (55) 



( Pringlelll98in . where ' = d/dr and the angular velocity is given in equation (|24|l . In the previous section, we found the tidal 
torque on the disc. However, this is only valid up to the position in the disc where the orbits cross. We adopt a torque function 
per unit disc mass on the disc of the form 

-Tz{r)(^Y, (56) 



dM ^ 

where g is & constant. To model the tidal effects of orbit crossings, we select the torque parameters so that the torque acts 
in the region where the particle orbits begin to cross. We choose (? = 4 so that the disc is truncated quickly. The form of Ta 
controls where the disc is truncated. We choose 

^ [0.5.^^7^ .^0.4.H 

r<0.4rH. 

To model the effects of the inflowing circumstellar gas, we model gas injection at some radius rinj over a narrow region of 
radial width 2w at a steady rate Afinj with the local Keplerian speed. Function Sir) describes the mass injection that we take 
to be 

Sir) = M., Hiir~r,^,)M^ ^^^^ 
27rrinj 2w 

where Hix) is unity for |a::| < 1 and zero otherwise. We adopt a width w = 0.0046 (rinj/rn)"^ r^. 



We solve this equation numerically on a fixed mesh that is uniform r 2 with 200 grid points (like iMartin et al.l 120071 ) 



We choose zero torque boundary conditions at the inner boundary nn — 10""^ rn and the outer boundary rout = 0.9 rn. The 
inner boundary allows the material there to be accreted by the planet. Its position has been chosen so that it is about equal 
to the radius of a Jupiter planet at a distance of 5 AU from a central solar mass star. With the tidal torque acting, the outer 
boundary is far enough out that it does not affect the evolution because the tidal torque prevents the mass from reaching 
the outer regions. However, if there is no tidal torque, the disc extends as far as it can and the outer boundary will affect 
the mass in the disc because mass is removed there. The viscous torque is given in equation 1)53^ and we take E = at the 
boundaries to have zero torques there. We initially take the surface density to be a constant but very small value and allow 
it to build up by mass accretion at a radius rjnj. 

We set the accretion rate on to the disc to be M = 1.78 x 10~* ilp. For a Jupiter mass planet orbiting a solar mass 
star at a radius of 5 AU this corresponds an accretion rate of M = 10~* Mq yr~^. We start with a disc of nearly zero mass. 
We note that changing the accretion rate does not change the results in Fig. [3l only the amount of mass in the disc. We 
parametrise the viscosity with the Q-prescription in equation ((5]). For a Keplerian disc with a — 10~^ and a constant disc 
aspect ratio H /r — 0.3, we have gas kinematic turbulent viscosity 

i/ = 1.56 X 10"'' rl^p. (59) 

Disc mass builds up as mass is injected. We ran the numerical code until it reached a steady state. In Fig. |3]we plot the 
surface density evolution of the disc for different injection radii. The surface density in each of the plots increases as time 
goes on. The top left plot shows the evolution without the tidal torque and an injection radius rinj — 0.2 rn. We see that the 
disc spreads out as far as it can out to the outer boundary, where mass is removed. All the other three plots in Fig. |3] include 
the tidal torque and we vary the position that mass is added from rinj = 0.1 rn (top right), 0.2 rn (bottom left) and 0.3 rn 
(bottom right). With the tidal torque acting, the disc cannot expand out to the outer boundary, it becomes truncated well 
inside that boundary. The position of the injection of the mass does not affect the outer boundary of the disc, it only mildly 
modifies the surface density profile of the steady state disc. 



6.5 One-Dimensional Analytic Solutions 

In this section we find stea dy-state analytical sol ution s for the surface density of a circumplanetary accretion disc. Our 
approach is similar to that of lCanup fc WardI (I2OO2I ) and I Ward fc Canup I J2OI0I ). but we include the strong tidal torques near 
the orbit crossing radius. There are three regions in the disc, as shown in Fig. ID Inside of the radius where mass is added, 
rinj, the disc acts as a normal accretion disc. Outside of this radius, the disc acts as a decretion disc in terms of the density 
profile, although gas does not actually fiow outward, as shown below. Instead, this region acts as mass reservoir. However, 
material in this region is subject to exchange with the interior region r < rinj by means of turbulent diffusion. A power law 
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Figure 3. Top left: The evolution of the disc without the tidal torque with rinj = 0.2 ru- The remaining plots show the evolution of 
the disc with the binary torque with rinj = 0.1 rpj (top right), ri^j = 0.2rH (bottom left) and rjnj = O.Srji (bottom right). The surface 
density increases in time. The times plotted, in order of increasing surface density in the plots, t = 0.01, 0.05, 0.18, 0.42, 0.67, 1.63, 5.50, 
20.96, 82.79 and 330.13 P where P is the orbital period of the planet, P = 27r/r2p. 



torque is applied only in the outer parts of the disc in r > rtrunc- We assume for simplicity that the disc is in exact Keplerian 
rotation about the planet with angular speed Q in the inertial frame. 
The accretion rate through the disc is 

M = 2nrur'£. (60) 
Notice that M is negative for accretion (ur < 0). Mass is injected at a radius r — rinj with rate Minj > 0, that is 

^ = M„j5(r - n„j) (61) 

and the angular momentum equation is 

^^^^^ = ^(T. + TgO + M„jr?„jni„j5(r - n„j), (62) 
where the viscous torque per unit radius (with equation I53p is given by 

= -STrr^nSi/ (63) 
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Figure 4. A sketch of the surface density about a planet in the three region model. In r < rjnj there is an accretion disc and in r > r^^j 
there is a decretion disc structure that acts as a mass reservoir. In this model, for r < rtrunc there is no imposed gravitational torque 
and for r > rtrunc a power law gravitational torque is applied. 



for a Keplerian angular velocity where fiinj = f2(rinj) and Tgr is the gravitational torque. Equation (|62|l is equivalent to 
equation (I55|l with the time derivative set to zero and injection width w set to (nearly) zero. 

In accordance with equation ((56]), the gravitational torque per unit radius on the disc is taken to be 



dTgr 

dr 



-27rr ET3(r) 



(64) 



where T3 is defined by equation (|57p . Integrating equation (|61ll we find 

-A'/inj + Cl 



M ■ 



(65) 



r < rinj 

Cl r > rinj, 

where Ci is a constant. With the boundary condition that M = at the disc outer edge, we see that Ci = 0. So there is no 
accretion flow outside of the injection radius, all of the material that is injected at rinj flows inwards to be accreted on to the 
planet. We now integrate equation (|62p to find 

C2 = 0, 
C2 = 0, 

C2 = 0, r> n 

where C2 is a constant to be determined by the inner boundary condition. As in Section [6.41 we have chosen the gravitational 
torque to be zero in r < rtmnc- 

In evaluating equation (|66p . we apply the form of kinematic viscosity given by equation (|59p . u{r) oc y'r, that holds for a 
constant a and constant H/r circumplanetary disc. With the zero viscous torque inner boundary condition at an inner radius 
Tin, equation (|66p requires C2 = — i\finjrfnr2in. We then obtain the standard viscous disc surface density in the inner parts of 



+ A/in 

+ Afinjrf„ 

Tv + Tgr + Minjn, 



f^inj 
^^inj 



r < rinj 
nj < r < rt rune 

rune ; 



(66) 



the disc r < rinj, designated by subscript 1, 



Si(r) 



Min 



1- — 



(67) 

we assume that 



I Pringldll981h . Therefore, Si(r) cx l/\/r, for Hnj > r S> Hn. In the other two regions of space where r > rinj, 
r S> rin, so that C2 is very small compared to the other torques. For simplicity, we take C2 = in these regions. 

In the region outside of where the mass is injected, but where no gravitational torque is applied, rinj < r < rtmn 
surface density is 

E2(r) = Ei(rinj)(^). 
In the outermost region, r ^ rtmnc, the surface density is 



E3 = E2 (rtrunc) ^— -) exp 



the 



(68) 



(69) 



3(1 4~ p) i^trunc^truncrtruncr _ 

where t'trunc = (^(rtmnc) and Sltrunc = n(rtrunc). We now have the surface density for the whole disc with the equations for 
El, E2 and E3. We compared these solutions to the one-dimensional cases of Section [6.41 with the same values for the input 
parameters. These steady state analytic solutions agree well (better than 6%) with the numerical solutions obtained at late 
times, plotted as the highest lines in the plots in Fig. (3] 
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Figure 5. Surface density distributions for a disc with constant a and constant H/r as a function of radius for cases where the mass 
is injected at radius r = 0.5. The surface density is normalized by its value at the injection radius. The dashed line is for a model that 
satisfies outer boundary condition T^/ = at the disc outer radius r^. The unit of radius on the plot is r^j in the case. The solid line 
is for a case that models the effects of the tidal field. The tidal field acts to sharply truncate the disc outside of the truncation radius 
^trunc- The density distribution in this case satisfies satisfies boundary condition I I70I1 and the unit of radius in the plot is rtrunc- The 
region r > rtrunc where the density drops to zero, is not shown. There are additional corrections to this profile due to the effects of gas 
pressure in the outer parts of the disc. 



Canup fc Ward 1 20021 ') have also considered the flow outside some radius where mass is injected, but did not include the 



effects of the tidal field. They have a normal accretion disc in the inner parts, as we have. But in the outer parts of the disc in 
their model, the material flows outwards and is removed from the disc at a rate 0.8 A/inj yHnj/'^ where is the outer edge of 
disc. For their standard model with nnj = 30rj and = 150rj (for Jupiter radius rj), we find the mass loss rate is 0.36 Minj. 
If we take the disc radius to be the orbit crossing radius, the mass loss rate is about 0.25 Afinj. They suggest that the removal 
of matter at the outer edge of the disc occurs through solar torques or collisions with highly shocked regions. However, in 
our model the tidal torque dominates the viscous torques in the outermost regions and prevents the outward flow of matter. 
Hence, in our picture of a steady state, the rate that mass is accreted on to the central planet is the same rate that it is 
injected into the disc. In addition, there is a difference in the density profiles between the models outside the mass injection 
radius. In that region, ICanup fc WardI l|2002l ) determined a gradual tapering that varies in radius as ^r^jr — 1 for disk radius 
rd. In our model, the density falls off more slowly in radius in the main body of the disc and is abruptly ter minated near the 



orbit crossing radius. The difference in the density structures can be traced to outer boundary conditions. I Canup fc Ward 



1 2OO2I) employed a zero viscous stress outer boundary condition, — 0. In our model, the effects of the tidal field can be 



approximately represented by an outer boundary condition at the tidal truncation radius (orbit crossing radius) that is given 
by 

T, = -Afinjrf^jfiinj, (70) 

as follows from the second relation of equation (|66p with C2 = 0. Fig. [5] comp ares the density distrib utions resulting from the 
two boundary conditions. The density distribution we obtain approaches the ICanup fc WardI (I2OO2I ) distribution in the limit 
that the injection radius goes to zero. There are further modifications to the density profile due to pressure effects that are 
described in Section [S] 

From equation (|69p . the width of the disc outer edge (ignoring pressure effects) is estimated as 

w = l.5unr ^ , (71) 

where the right-hand side is evaluated near the outer edge. If the disc flow is smooth and there are no resonances, the phase 
lag that produces a tidal torque due to the turbulent viscosity is estimated as 86 ~ i//(r^r2), the inverse Reynolds number of 
the flow. The resulting torque per unit disc mass is then estimated as \dT^^ j dM\ ~ |<[>2p/(r^f2^), where $2 is the m = 2 
component of the tidal potential. Using equation (|71() and assuming |$2| ^ r^V^ , it follows that w 3> r. Therefore, this torque 
is not capable of truncating the disc. Instead, the much stronger tidal torque, due to orbit crossings and certain resonance 
effects discussed in Section [S] truncates the disc in a small region of space. Torques from the region well inside of where orbit 
crossings occur are then too weak to truncate the disc. 

The analytic solutions can also be considered to apply when the matter being injected into the disc has a nonKeplerian 
velocity, that is fiinj at rinj is nonKeplerian. The gas in such a case adds angular momentum per unit mass to the system at a 
rate that differs from the Keplerian rate. In reality, such gas would undergo a strong Kelvin Helzholtz instability with the disc. 
Since the gas being added to the disc has a much lower density than the disc, we would expect this gas to be rapidly entrained 
by the disc that remains in nearly Keplerian rotation. An example of such a situation occurs when coplanar inflowing gas 
meets the circumplanetary disc at the disc outer edge. We expect the inflowing gas in this case to be subKeplerian. In such 
a case, the decretion disc or mass reservoir region, as depicted in Fig. [4] is very small. Instead, nearly entire disc behaves an 
accretion disc. The subKeperian injected gas can be considered to impart a negative torque on an otherwise Keplerian disc, 
in addition to the tidal torque. On the other hand, if much of the inflowing gas flows over the disc before becoming entrained 
by the disc, then the injection radius will be smaller than the disc outer radius, as in the cases plotted in Fig. [S] 
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Figure 6. Average specific angular momentum (angular momentum per unit mass) (J) in a circumplanetary disc (total disc angular 
momentum divided by disc mass), based on the steady state model of Section 6.5, plotted as a function of the specific angular momentum 
of the accreting gas Jinj- The quantities on the horizontal and vertical axes are normalized by r^Qp. The range of the horizontal axis 
extends to Ji^j at the disc outer edge, r = 0.41rH. The results show that the average specific disc angular momentum is nearly independent 
of the specific angular momentum of the accreting gas, even for small values. This model assumes that the gas fiows over the disc and 
is entrained vertically at the local Keplerian speed for the mass injection radius. If planar subKeperian gas is added at the disc outer 
edge, we expect the variations of (J) to be even smaller, since the gas distribution within the main body of the disk behaves as a simple 
accretion disc. 




Beyond the details of the flow properties in particular models, the results suggest that the disc structure is generally 
largely independent of the angular momentum per unit mass of the inflowing gas, as seen in Figs. [3] and [6] The disc density 
distribution is mildly influenced by this quantity. The radial derivative of the surface density undergoes order unity changes 
at rinj. In the model considered this section, the logarithmic surface density gradient changes from -1/2 to -1 across nnj. Also, 
for a flxed mass injection rate, the density in the outer parts of the disc scales as \/ nnj/^trunc- But such variations do not 
lead to strong changes in disc structure by typically plausible values of rinj > 0.05rH. This is particularly true if the injection 
radius is near the disc outer edge. Instead, the disc structure interior to the disc tidal truncation region is dominated by the 
effects of disc turbulent viscosity. 



7 DISC MODEL WITH PRESSURE 



To include the effects of gas pressure, we modeled the disc using a two-dimensional SPH code with 10^ particles. The SPH 
parameter asPH was set to unity in these simulations, while /3sph was set to zero. The disc sound speed was crudely modeled as 
a constant equal to 0.3 /i^^^aQp = 0.43 ilpru throughout. The disc aspect ratio was H/r ~ 0.16 at i? = OARh (in dimensionless 
units of these equations Ru = 3"^''^), where free particle orbit crossings occur. The fluid equations were taken in the Hill 
approximation, following the force equation (|30|) . The particles were initially distributed so that surface density was initially 
E(r, 6) oc 1/R in an annulus Rin < R < i?out- We chose Rin = 0.07i?H and have two models, one with i?out = 0.37iH and the 
second with -Rout = 0.6i?H which we plot in Fig. [T] The initial velocities were taken to be circular rotation having angular 
speed f2 — Sip in the corotating frame, with Q, given by equation ((24]). Particles are removed from the simulation if they reach 
the inner boundary at R = 0.07i?H or the outer boundary located at J? = 0.86J?h. For each particle removed, a particle is 
injected at a random angle and random radius between R = 0.22i?H and R = 0.36i?H. Therefore, the number of particles is 
flxed at 10^ at all times. 

Fig. [8] plots the angular momentum evolution as a function of time for the two initial disc sizes. The smaller disc starts 
with much less angular momentum than the larger one. But after only about 3 planetary orbits both discs have similar values 
of angular momentum and approach a steady state. In Fig. |9]we plot the particles at a time of 6 planet orbital periods. We 
see that particle distributions look very similar and the disc has reached a nea r steady state. The viscous ti mescale of the disc 



can be estimated as ~ r /v. In SPH, we have that a ~ O.laspH = 0.1 (e.g., lArtvmowicz fc Lubowlll994l ). The value of the 



kinematic viscosity v is estimated as 0.1aspn{H/r)^r^fl and the viscous timescale evaluates to about 8 orbits at i? ~ 0.3i?H- 
The simulated disc has therefore settled to a near equilibrium state on a timescale that is of order the estimated viscous 
timescale. This ki nematic viscosity is a factor of ten or more larger than what is typically taken taken in circumplanetary disc 
simulations (e.g.. IP'Angelo et al.ll20oi : lAvliffe fc Bat j|2009l ). Consequently, in those simulations, the timescales to reach a 
steady state are longer than the case here by similar factor. (In addition, previous simulations have generally not start e d wit h 
a circumplanetary disc. Some time is required for its formation from inflowing gas.) For the simulations in iMachidal (|2009l ). 
the timescale would be determined by the inherent viscosity in the code due to the flnite differencing, since no viscous terms 
were included. 

The phasing of the arms determines the sign of the gravitational torque in the different regions of space as shown in 
Fig. 1101 The left plot shows the cumulative gravitational torque as a function of radius in a steady state disc after 6 planetary 
orbits. On the right is a plot of the shape of the spiral in the disc. The sign of the gravitational torque is equal to sign of 
~xy and so is negative in the flrst and third quadrants the torque and positive in the second and fourth. As seen in the left 
plot, the gravitational torque is negative in the outermost region of the disc where it is strongest. The spiral there lies in the 
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Figure 7. The initial configuration for tiie two SPH runs. Left: the particles are randomly distributed in radius of O.OTrfj < r < 0.3rH. 
Right: The particles are randomly distributed randomly in angle and randomly in radius within and annulus defined by COTrj] < r < 
0.65rH. 
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Figure 8. The angular momentum evolution in time for the initially small disc (left) and the large disc (right), plotted in Fig. [7] The 
units of angular momentum are M^r^P~^ , for a disc of mass M^. The unit of time is the planet orbit period, P. 



first and third quadrants, as seen in Fig. [9] in order the provide a net negative gravitational torque that removes angular 
momentum from the disc as gas accretes. 

The left plot in Fig.[TO]aIso shows the cumulative disc mass distribution. More than 90% of the disc mass is located inside 
a radius of 0.41rH, where orbits crossings of free particles occur. The local gravitational torque is negative for r > 0.35rH. On 
the other hand, the cumulative gravitational torque on the disc inside the orbit crossing radius is actually positive. This does 
not mean that the gas gains angular momentum in this region. The viscous torques compensate by transporting the angular 
momentum out of this region so that the overall torque is negative. The negative cumulative gravitational disc torque is 
achieved somewhat outside the orbit crossing radius. A small amount of disc mass, less than 10%, that is located outside that 
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Figure 9. The particle distribution in space for the two SPH simulations corresponding to those with the initial configuration in Fig. [7] 
at a time of 6 planetary orbital periods. 
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Figure 10. Left: Cumulative disc mass (solid) in units of total disc mass Md and cumulative gravitational torque (dashed) on the disc 
in units of 0.01M(jr^f2p starting at the disc center as a function of radius in units of rji at a time of 6 planetary orbits. Right: Plot of the 
particle positions in units of rji with the point size weighted by density, in order to enhance the visibility of the spiral wave in the main 
body of the disc. The sign of the gravitational torque in each quadrant is also shown. At small radii, r ~ 0.2, the spiral lies in the first 
and third quadrants where the gravitational torque on the disc is negative, and the cumulative gravitational torque begins as negative 
in the left panel. At r ~ 0.3, the spiral lies in the second and fourth quadrants where the gravitational torque on the disc is positive, 
and the cumulative gravitational torque turns positive. The outer parts of the spiral lie in the first and third quadrants and provide a 
negative gravitational torque contribution that causes the total torque on the disc to be negative. 



radius provides the net negative torque on the disc. In that region of space, r > 0.41rH, the spiral arms are very prominent in 
Figure[9]and gas response is hi ghly nonlinear. A sim ilar spiral structure was found in the case of a warm disc with H/r = 0.1 
within a binary star system bv lSavoniie et al.l (|l994l ) that we have argued should be similar to the circumplanetary disc case 
(see Section |4]). 

The SPH code is less accurate in these outer low density regions because the interparticle pressure force calculation is 
more approximate. Also, we have simplified the infiow on to the disc to be occurring within an annulus in the disc. But, 
we expect that the general properties of gravitational torque to generally hold. That is, the gravitational torque involves a 
relatively small amount of gas in the outer parts of disc. 



8 DISCUSSION AND CONCLUSIONS 



We have analysed the dynamics of a circumplanetary disc as an accretion disc subject to the tidal forces of the central star. 
We applied several techniques: ballistic particle orbits, one-dimensional simulations and analytic models, and two-dimensional 
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Figure 11. Angular momentum per unit disc mass, J, in units of f^pf^ plotted against distance from the planet in units of rjj on a 
log-log scale. The solid line is obtained from the SPH simulations described in Section 7 at a time of 6 orbits, while the dashed line is 
the Keplerian case. The solid curve is based on binning in 30 radial bins, in order to smooth fluctuations. For r ~ 0.3, inside of where 
orbits cross (r ~ 0.4), J in the disc departs from the Keplerian case due to the effects of gas pressure as the disc tapers. The outer disc 
velocity is somewhat subKeplerian because the pressure forces act radially outward. 




SPH simulations. We have shown that the gas dynamical equations for a circumplanetary disc can be rescaled to a form that 
is similar to those previously analysed for a disc in a binary star system with mass ratio of unity, but with a modified (Hill) 
potential. The size of a circumplanetary disc is determined by the requirement that tidal torques from the central star remove 
angular momentum at a rate required by steady inflow on to the planet. In the crudest approximation that is accurate for a 
cold disc, disc streamlines consist of a set of stable, nested periodic ballistic particle orbits (see Fig. [!}. For a cold disk, the 
disc radius is determined to be the location of orbit intersections or instability that occurs at 0.41rH, for Hill radius rn. This 
radius is close to the effective disc radius (inside of which ~ 90% of the mass is located) found in the SPH simulations for a 
warm disc (see Figs. [9l and 1101) . In principle, the observational determination of the size of this disc then provides a constraint 
on the planet mass. 

In some previous disc models, the angular momentum of the inflowing gas was thought to control the disc size. According 
to these models, the compactness of the giant planet satellite systems, within a region of radius much less than ru, is a 
consequen ce of a low angular momentum inflow that penetrates well inside the Hill sphere and depo sits presatellite material 
there (e.g. lLissauer1ll995l : [Canup fc Wardlbood : [Mosaueira fc Estradal [200^ : IWard fc Canup Iboiol ). We suggest a different 
physical model for the processes that control the circumplanetary disc size and structure. We find that although the angular 
momentum per unit disc mass (specific angular momentum) of the inflowing gas plays a role in determining whether a 
circumplanetary disc can form, it has a minor influence on the disc structure. A steady state accretion disc does not achieve 
as compact a form as desired by the se pr evious models. These results imply that the common scale of the disc seen in the 
simulations bv lAvliffe fc Bate! (|2009l ) and lD'Angelo et alj I (|2002h is not due to the specific angular momentum of the inflow, 
but instead to the tidal truncation effects. The angular momentum within the disc is redistributed by viscous torques. The 
disc structure inside the truncation radius is mainly determined by the properties of the disc turbulent viscosity, rather than 
the specific angular momentum of the inflowing gas (see uppermost curves in Fig. [S]). The average specific angular momentum 
in the disc is nearly independent of the speciflc angular momentum of the accreting gas (see Fig. (Gjl . Tidal effects have a 
large-scale influence on the disc density distribution, as discussed in Section [6.51 

We have concentrated on the case that the planet opens a gap in the disc, as is expected for a Jupiter mass planet. Some 
previous studies sug gested that a small disc of size ~ rH/48 could result prior to and during early stages of gap opening (e.g. 
Estrada et al. 20081 ) . Our results suggest that a turbulent accretion disc would be truncated at the orbit crossing radius, even 



if the inflowing gas ha s low specific angular momentum. There are are few qualiflcations in this statement. The simulations of 
Avliffe fc Ba"t3 l|2009l ) found that circumplanetary discs do not exist prior to gap opening (planet masses less than 100 Earth 
masses). If low angular momentum gas is accreted in early stages of gap opening, then pressure forces from within the Hill 
sphere could prevent a disc from extending out to the orbit crossing radius, if the planet's Bondi radius is sufficiently smaller 
than its Hill radius. Under such conditions, we would expect the disc to expand to the Bondi radius by viscous torques, since 
tidal torques would not be adequate to truncate the disc inside the Bondi radius. The disc would then need to lose its angular 
momentum to the surrounding gas in the circumstellar disc. Whether such a configuration is possible is unclear. In addition, 
gap opening or planet envelope contraction (below its Hill sphere) might not begin when the Bondi radius is smaller than the 
Hill radius and so a disc may not form under such conditions. 

During the T Tauri accretion phase, circumplanetary disc aspect ratios are expected to be large, H/r ~ 0.3, as a 
consequence of the weak gravitational forces due to the planet and relatively high temperatures (see equation Q). There are 
some implications of these high disc aspect ratios. One is that the disc edge tapering is more gradual than in the thin disc 
case. The tapering must be on a scale > H in order to avoid instabilities that occur for sharper edge s (e.g. Yang fc Menoul 
The tapering can then have an influence on a substantial portion of the circumplanetary disc. lAvliffe fcBatX([200^ 



estimated the disc outer radius based on where the angular momentum per unit disc mass departs from Keplerian and begins 
to decline in radius. They obtained a disc outer radius of ~ 0.35rH (see the upper right panel of their Fig. 2). This radius 
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is somewhat inside the orbit crossing radius of ~ 0.41r-H. The effect of the disc tapering on the distribution of disc angular 
momentum in our SPH simulations can be see n in Fig. 1111 It leads to subKeplerian velocities inside the orbit crossing radius, 
similar to the results of Avliffe fc Batel ( 20091 '). The results of these two studies can be further reconciled by considering the 
differences in the models. In their simulations, the disc is hotter than in our simulations. Consequently the pressure forces 
and their effects on departures from Keplerian rota tion in the out e r par ts of the disc should be stronger than in Fig. Illl 
as is consistent with their results. Also the work of I Avliffe fc Bate l|2009l ) includes the inflowing gas from larger radii that 
is omitted in our simulations. This inflowing subKeplerian gas and its interaction with the disc outer edge as it becomes 
entrained could further modify the results plotted in Fig. Illl It is also possible that the disc radius is affected by the location 
of wave damping, as discussed below. Such effects should be explored further. 

Another implication of the high values of H/r is that a small fraction of the disc mass can reside at larger disc radii than 
the orbit-crossing radius and so be subject to the effects of resonances that lie outside this radius, as discussed in Section [S] 
(see Fig. llOf) . In addition, off-resonant tidal forcing of the outer parts of the disc that lie inside the orbit crossing radius 
can play an important role. The resonance width depends on {H/r)'^^^ and is not small. The substantial width can allow 
the resonant region to overlap with the denser parts of disc. As a result of these two effects, two-armed spiral waves can be 
launched in such discs, even though exact resonances do not lie within the main body of the discs. We have found evidence 
of such waves in our SPH simulations (see Fig. [9]). For these waves to play a role in extracting angular momentum from the 
disc, they must damp in order to introduce irreversibility. Otherwise, a standing wave is produced that results in little or no 
torque on the disc. Some damping can be produced by the disc turbulent viscosity in the fractional amount that is roughly 
the ratio of the viscous wave damping rate to the wave group propagation rate, about a{r/H). This ratio can be quite small 
for a <^ H/r and the torque is reduced to a similarly small fraction of its potentially maximum value. For a thin disc, this 
damping could occur as launched waves propagate, steepen and shock, since their wavelengths are short compared to the disc 
radius. But for warm circumplanetary discs, this process is less important because the wavelengths are not short compared to 
the disc radius. In the SPH simulations of Section 7, the damping seems to occur from the strongly nonlinea r forcing in the 
outer disc, similar to the case of mildly warm discs {H/r = 0.1) previously investigated for binary star systems (jSavoniie et al 



1994 ). The SPH simulations suggest that the negative torque is produced in the outer parts of the disc, somewhat inside and 



beyond the radius where nested particle orbits cross or become unstable. This torque involves a relatively small amount of 
disc mass (see Fig. I10|l. However, the discs in our SPH simulations were two-dimensional and not as warm as could occur, 
where H/r ~ 0.3. The wave damping issue should be explored in future work. 

The accretion disc model does not provide an obvious explanation for the locations of the regular satellites of Jupiter and 
Saturn that occur within 0.06rH of the planet. As discussed above, the disc structure is insensitive to the angular momentum 
of the inflowing gas and satellites lie well inside the tidal truncatio n radiu s of the disc. This conclus ion is consistent with 
the results of recent three-dimensional simulations of lEstrada et al.1 (|2008l l and lAvliffe fc Batel (|2009f ). We have considered 
circumplanetary discs with a smoothly varying turbulent viscosity with radius. But such discs could also harbor dead zones, 
as discussed in Section (2] Such situations could result in rapid disc density variations at their boundaries. The existence of 
satellites requires survival against the effects of migration. Such density varia tions could affect satellite migration and possibly 
trapping satellites, since migration rates depend on density gradients (e.g., Matsumura. Pudritz fc Thomme^ 2007 1. But, if 
the inner dead zone boundary is due to thermal ionisation at a temperature of 10^ K, then the disc temperatures there would 
be too hot for the survival of an icy satellite. Of course, other sources of ionisation could change the locations of dead zone 
edges. In addi tion, the nonlinear feedback from waves in a low viscosity disc could slow migration jWardll 19971 : lRafikovll2o"o3l : 
Li et al. 20091 ) . However, this slowing is less effective for the warmer circumplanetary discs. If the satellite formation occurs 



after partial disc depletion, as suggested by Canup fc Ward (2002), then the ionisation through the vertical extent of the disc 
becomes easier and the dead zones are less likely to occur at that stage. 

Circumplanetary discs are not very bright as we see in equation ^ . Time dependent accretion and outbursts could occur 
in de ad zones, as has been suggested to explain FU Orionis outbursts in young stellar systems (jArmitage. Livio fc Pringle 



20011 ). During a circumplanetary disc outburst, the ratio of the circumplanetary to circumstellar luminosities could be much 
higher. 
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